//||=========================================================================||
//||                                                                         ||
//||                                main.cpp                                 ||
//||                                                                         ||
//||                           by the_best_goose))                           ||
//||                                                                         ||
//||=========================================================================||

#include <stdio.h>
#include <assert.h>
#include "MolDynSys.h"

//-----------------------------------------------------------------------------

int main()
{
    const int particle_count = 512;
    const int steps_count = 1000;
    MolDynSys sys(particle_count, 0.8, NO_TEST);
    sys.InitSystem();

    for(int i = 0; i < 10; ++i)
    {
        sys.RunThermostat(1.24, 0.001, 0.01);
        sys.Relax(0.02, 0.01);
    }
    float energy = sys.GetPotEnergy();
    printf("Default potential energy per particle = %f\n", energy / particle_count);


    float step = 0.01;
    double average_energy = 0.0;
    double average_kin_energy = 0.0;
    
    for(int i = 0; i < steps_count; ++i)
    {
        sys.NextStep(step);
        sys.Print();

        printf("t = %f\r", i * step);
        average_energy += sys.GetPotEnergy();
        average_kin_energy += sys.GetKinEnergy();
        fflush(stdout);
    }
    average_energy /= particle_count * steps_count;
    average_kin_energy /= particle_count * steps_count;
    
    printf("Average Potential energy per particle = %f\n", average_energy);
    printf("Average Temperature = %f", average_kin_energy / 1.5);

    return 0;
}
